home *** CD-ROM | disk | FTP | other *** search
/ Languguage OS 2 / Languguage OS II Version 10-94 (Knowledge Media)(1994).ISO / gnu / glibc108.zip / glibc108 / stdio / printf_fp.c < prev    next >
C/C++ Source or Header  |  1994-05-10  |  18KB  |  737 lines
  1. /* Floating-point printing for `printf'.
  2.    This is an implementation of a restricted form of the `Dragon4'
  3.    algorithm described in "How to Print Floating-Point Numbers Accurately",
  4.    by Guy L. Steele, Jr. and Jon L. White, presented at the ACM SIGPLAN '90
  5.    Conference on Programming Language Design and Implementation.
  6.  
  7. Copyright (C) 1992, 1993, 1994 Free Software Foundation, Inc.
  8. This file is part of the GNU C Library.
  9.  
  10. The GNU C Library is free software; you can redistribute it and/or
  11. modify it under the terms of the GNU Library General Public License as
  12. published by the Free Software Foundation; either version 2 of the
  13. License, or (at your option) any later version.
  14.  
  15. The GNU C Library is distributed in the hope that it will be useful,
  16. but WITHOUT ANY WARRANTY; without even the implied warranty of
  17. MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  18. Library General Public License for more details.
  19.  
  20. You should have received a copy of the GNU Library General Public
  21. License along with the GNU C Library; see the file COPYING.LIB.  If
  22. not, write to the Free Software Foundation, Inc., 675 Mass Ave,
  23. Cambridge, MA 02139, USA.  */
  24.  
  25. #include <ansidecl.h>
  26. #include <ctype.h>
  27. #include <stdio.h>
  28. #include <float.h>
  29. #include <limits.h>
  30. #include <math.h>
  31. #include <stdarg.h>
  32. #include <stdlib.h>
  33. #include <localeinfo.h>
  34.  
  35. #include <printf.h>
  36.  
  37. #define NDEBUG
  38. #include <assert.h>
  39.  
  40. #define    outchar(x)                                  \
  41.   do                                          \
  42.     {                                          \
  43.       register CONST int outc = (x);                          \
  44.       if (putc (outc, s) == EOF)                          \
  45.     return -1;                                  \
  46.       else                                      \
  47.     ++done;                                      \
  48.     } while (0)
  49.  
  50. #if FLT_RADIX != 2
  51.  #error "FLT_RADIX != 2.  Write your own __printf_fp."
  52. #endif
  53.  
  54. #undef alloca            /* gmp-impl.h defines it again.  */
  55. #include "gmp.h"
  56. #include "gmp-impl.h"
  57. #include "longlong.h"
  58.  
  59. #ifndef NDEBUG
  60. static void mpn_dump (const char *str, mp_limb *p, mp_size_t size);
  61. #define MPN_DUMP(x,y,z) mpn_dump(x,y,z)
  62. #else
  63. #define MPN_DUMP(x,y,z)
  64. #endif
  65.  
  66. extern mp_size_t __mpn_extract_double (mp_ptr res_ptr, mp_size_t size,
  67.                        int *expt, int *is_neg,
  68.                        double value);
  69.  
  70. /* 1000 words is an arbitrarily chosen size that is big enough.  Some day I
  71.    should figure out how big it actually needs to be (which should be
  72.    computable from the <float.h> parameters).  The size initialization is
  73.    just a notable random number to catch uninitialized variable bugs.  */
  74. #define MPN_VAR(name) mp_limb name[1000]; mp_size_t name##size /*= 314159265*/
  75. #define MPN_ASSIGN(dst,src) \
  76.   memcpy (dst, src, (dst##size = src##size) * sizeof (mp_limb))
  77. #define MPN_POW2(dst, power) \
  78.   do {                                          \
  79.     MPN_ZERO (dst, (power) / BITS_PER_MP_LIMB);                      \
  80.     dst[(power) / BITS_PER_MP_LIMB] =                          \
  81.       (mp_limb) 1 << (power) % BITS_PER_MP_LIMB;                  \
  82.     dst##size = (power) / BITS_PER_MP_LIMB + 1;                      \
  83.   } while (0)
  84.  
  85. /* Compare *normalized* mpn vars.  */
  86. #define MPN_GT(u,v) \
  87.   (u##size > v##size || (u##size == v##size && __mpn_cmp (u, v, u##size) > 0))
  88. #define MPN_LT(u,v) \
  89.   (u##size < v##size || (u##size == v##size && __mpn_cmp (u, v, u##size) < 0))
  90. #define MPN_GE(u,v) \
  91.   (u##size > v##size || (u##size == v##size && __mpn_cmp (u, v, u##size) >= 0))
  92. #define MPN_LE(u,v) \
  93.   (u##size < v##size || (u##size == v##size && __mpn_cmp (u, v, u##size) <= 0))
  94. #define MPN_EQ(u,v) \
  95.   (u##size == v##size && __mpn_cmp (u, v, u##size) == 0)
  96. #define MPN_NE(u,v) \
  97.   (!MPN_EQ(u,v))
  98.  
  99. int
  100. DEFUN(__printf_fp, (s, info, args),
  101.       FILE *s AND CONST struct printf_info *info AND va_list *args)
  102. {
  103.   mp_limb cy;
  104.  
  105.   int done = 0;
  106.  
  107.   /* Decimal point character.  */
  108.   CONST char *CONST decimal = _numeric_info->decimal_point;
  109.  
  110.   LONG_DOUBLE fpnum;        /* Input.  */
  111.   int is_neg;
  112.  
  113.   MPN_VAR (f);            /* Fraction.  */
  114.  
  115.   int e;            /* Base-2 exponent of the input.  */
  116.   CONST int p = DBL_MANT_DIG;    /* Internal precision.  */
  117.   MPN_VAR (scale); MPN_VAR (scale2); MPN_VAR (scale10); /* Scale factor.  */
  118.   MPN_VAR (loerr); MPN_VAR (hierr); /* Potential error in the fraction.  */
  119.   int k;            /* Digits to the left of the decimal point.  */
  120.   int cutoff;            /* Where to stop generating digits.  */
  121.   MPN_VAR (r); MPN_VAR (r2); MPN_VAR (r10); /* Remainder.  */
  122.   int roundup;
  123.   int low, high;
  124.   char digit;
  125.  
  126.   MPN_VAR (tmp);        /* Scratch space.  */
  127.  
  128.   int j;
  129.  
  130.   char type = tolower (info->spec);
  131.   int prec = info->prec;
  132.   int width = info->width;
  133.  
  134.   /* This algorithm has the nice property of not needing a buffer.
  135.      However, to get the padding right for %g format, we need to know
  136.      the length of the number before printing it.  */
  137.  
  138. #ifndef    LDBL_DIG
  139. #define    LDBL_DIG    DBL_DIG
  140. #endif
  141. #ifndef    LDBL_MAX_10_EXP
  142. #define    LDBL_MAX_10_EXP    DBL_MAX_10_EXP
  143. #endif
  144.  
  145.   char *buf = __alloca ((prec > LDBL_DIG ? prec : LDBL_DIG) +
  146.             LDBL_MAX_10_EXP + 3); /* Dot, e, exp. sign.  */
  147.   register char *bp = buf;
  148. #define    put(c)    *bp++ = (c)
  149.  
  150.  
  151.   /* Produce the next digit in DIGIT.
  152.      Return nonzero if it is the last.  */
  153.   inline int hack_digit (void)
  154.     {
  155.       int cnt;
  156.       mp_limb high_qlimb;
  157.  
  158.       --k;
  159.       cy = __mpn_mul_1 (r10, r, rsize, 10);
  160.       r10size = rsize;
  161.       if (cy != 0)
  162.     r10[r10size++] = cy;
  163.  
  164.       MPN_DUMP ("r", r, rsize);
  165.       MPN_DUMP ("r10", r10, r10size);
  166.       MPN_DUMP ("scale", scale, scalesize);
  167.  
  168.       /* Compute tmp = R10 / scale  and  R10 = R10 % scale.  */
  169.       count_leading_zeros (cnt, scale[scalesize - 1]);
  170.       if (cnt != 0)
  171.     {
  172.       mp_limb norm_scale[scalesize];
  173.       mp_limb cy;
  174.       assert (scalesize != 0);
  175.       __mpn_lshift (norm_scale, scale, scalesize, cnt);
  176.       assert (r10size != 0);
  177.       cy = __mpn_lshift (r10, r10, r10size, cnt);
  178.       if (cy != 0)
  179.         r10[r10size++] = cy;
  180.       high_qlimb = __mpn_divmod (tmp, r10, r10size, norm_scale, scalesize);
  181.       tmp[r10size - scalesize] = high_qlimb;
  182.       r10size = scalesize;
  183.       __mpn_rshift (r10, r10, r10size, cnt);
  184.     }
  185.       else
  186.     {
  187.       high_qlimb = __mpn_divmod (tmp, r10, r10size, scale, scalesize);
  188.       tmp[r10size - scalesize] = high_qlimb;
  189.       r10size = scalesize;
  190.     }
  191.  
  192.       MPN_DUMP ("high_qlimb", &high_qlimb, 1);
  193.       MPN_DUMP ("r10", r10, r10size);
  194.  
  195.       /* We should have a quotient < 10.  It might be stored */
  196.       high_qlimb = tmp[0];
  197.       digit = '0' + high_qlimb;
  198.  
  199.       r10size = __mpn_normal_size (r10, r10size);
  200.       if (r10size == 0)
  201.     /* We are not prepared for an mpn variable with zero limbs.  */
  202.     r10size = 1;
  203.  
  204.       MPN_ASSIGN (r, r10);
  205.       assert (rsize != 0);
  206.       cy = __mpn_lshift (r2, r, rsize, 1);
  207.       r2size = rsize;
  208.       if (cy != 0)
  209.     r2[r2size++] = cy;
  210.  
  211.       cy = __mpn_mul_1 (loerr, loerr, loerrsize, 10);
  212.       if (cy)
  213.     loerr[loerrsize++] = cy;
  214.       cy = __mpn_mul_1 (hierr, hierr, hierrsize, 10);
  215.       if (cy)
  216.     hierr[hierrsize++] = cy;
  217.  
  218.       low = MPN_LT (r2, loerr);
  219.  
  220.       /* tmp = scale2 - hierr; */
  221.       if (scale2size < hierrsize)
  222.     high = 1;
  223.       else
  224.     {
  225.       cy = __mpn_sub (tmp, scale2, scale2size, hierr, hierrsize);
  226.       tmpsize = scale2size;
  227.       high = cy || (roundup ? MPN_GE (r2, tmp) : MPN_GT (r2, tmp));
  228.     }
  229.  
  230.       if (low || high || k == cutoff)
  231.     {
  232.       /* This is confusing, since the text and the code in Steele's and
  233.          White's paper are contradictory.  Problem numbers:
  234.          printf("%20.15e\n", <1/2^106>) is printed as
  235.          1.232595164407830e-32 (instead of 1.232595164407831e-32)
  236.          if we obey the description in the text;
  237.          1/2^330 is badly misprinted if we obey the code.  */
  238.       if (high && !low)
  239.         ++digit;
  240. #define OBEY_TEXT 1
  241. #if OBEY_TEXT
  242.       else if (high && low && MPN_GT (r2, scale))
  243. #else
  244.       else if (high == low && MPN_GT (r2, scale))
  245. #endif
  246.         ++digit;
  247.       return 1;
  248.     }
  249.  
  250.       return 0;
  251.     }
  252.  
  253.   const char *special = NULL;    /* "NaN" or "Inf" for the special cases.  */
  254.  
  255.   /* Fetch the argument value.  */
  256.   if (info->is_long_double)
  257.     fpnum = va_arg (*args, LONG_DOUBLE);
  258.   else
  259.     fpnum = (LONG_DOUBLE) va_arg (*args, double);
  260.  
  261.   /* Check for special values: not a number or infinity.  */
  262.  
  263.   if (__isnan ((double) fpnum))
  264.     {
  265.       special = "NaN";
  266.       is_neg = 0;
  267.     }
  268.   else if (__isinf ((double) fpnum))
  269.     {
  270.       special = "Inf";
  271.       is_neg = fpnum < 0;
  272.     }
  273.  
  274.   if (special)
  275.     {
  276.       int width = info->prec > info->width ? info->prec : info->width;
  277.  
  278.       if (is_neg || info->showsign || info->space)
  279.     --width;
  280.       width -= 3;
  281.  
  282.       if (!info->left)
  283.     while (width-- > 0)
  284.       outchar (' ');
  285.  
  286.       if (is_neg)
  287.     outchar ('-');
  288.       else if (info->showsign)
  289.     outchar ('+');
  290.       else if (info->space)
  291.     outchar (' ');
  292.  
  293.       {
  294.     register size_t len = 3;
  295.     while (len-- > 0)
  296.       outchar (*special++);
  297.       }
  298.  
  299.       if (info->left)
  300.     while (width-- > 0)
  301.       outchar (' ');
  302.  
  303.       return done;
  304.     }
  305.  
  306.   /* Split the number into a fraction and base-2 exponent.  The fractional
  307.      part is scaled by the highest possible number of significant bits of
  308.      fraction.  We represent the fractional part as a (very) large integer. */
  309.  
  310.   fsize = __mpn_extract_double (f, sizeof (f) / sizeof (f[0]),
  311.                 &e, &is_neg, fpnum);
  312.  
  313.   if (prec == -1)
  314.     prec = 6;
  315.   else if (prec == 0 && type == 'g')
  316.     prec = 1;
  317.  
  318.   if (type == 'g')
  319.     {
  320.       if (fpnum != 0)
  321.     {
  322.       if (is_neg)
  323.         fpnum = - fpnum;
  324.  
  325.       if (fpnum < 1e-4)
  326.         type = 'e';
  327.       else
  328.         {            /* XXX do this more efficiently */
  329.           /* Is (int) floor (log10 (FPNUM)) >= PREC?  */
  330.           LONG_DOUBLE power = 10;
  331.           j = prec;
  332.           if (j > p)
  333.         j = p;
  334.           while (--j > 0)
  335.         {
  336.           power *= 10;
  337.           if (fpnum < power)
  338.             /* log10 (POWER) == floor (log10 (FPNUM)) + 1
  339.                log10 (FPNUM) == Number of iterations minus one.  */
  340.             break;
  341.         }
  342.           if (j <= 0)
  343.         /* We got all the way through the loop and F (i.e., 10**J)
  344.            never reached FPNUM, so we want to use %e format.  */
  345.         type = 'e';
  346.         }
  347.     }
  348.  
  349.       /* For 'g'/'G' format, the precision specifies "significant digits",
  350.      not digits to come after the decimal point.  */
  351.       --prec;
  352.     }
  353.  
  354.   if (fsize == 1 && f[0] == 0)
  355.     /* Special case for zero.
  356.        The general algorithm does not work for zero.  */
  357.     {
  358.       put ('0');
  359.       if (tolower (info->spec) != 'g' || info->alt)
  360.     {
  361.       if (prec > 0 || info->alt)
  362.         put (*decimal);
  363.       while (prec-- > 0)
  364.         put ('0');
  365.     }
  366.       if (type == 'e')
  367.     {
  368.       put (info->spec);
  369.       put ('+');
  370.       put ('0');
  371.       put ('0');
  372.     }
  373.     }
  374.   else
  375.     {
  376.       cutoff = -prec;
  377.  
  378.       roundup = 0;
  379.  
  380.       if (e > p)
  381.     {
  382.       /* The exponent is bigger than the number of fractional digits.  */
  383.       MPN_ZERO (r, (e - p) / BITS_PER_MP_LIMB);
  384.       if ((e - p) % BITS_PER_MP_LIMB == 0)
  385.         {
  386.           MPN_COPY (r + (e - p) / BITS_PER_MP_LIMB, f, fsize);
  387.           rsize = fsize + (e - p) / BITS_PER_MP_LIMB;
  388.           assert (rsize != 0);
  389.         }
  390.       else
  391.         {
  392.           assert (fsize != 0);
  393.           cy = __mpn_lshift (r + (e - p) / BITS_PER_MP_LIMB, f, fsize,
  394.                  (e - p) % BITS_PER_MP_LIMB);
  395.           rsize = fsize + (e - p) / BITS_PER_MP_LIMB;
  396.           if (cy)
  397.         r[rsize++] = cy;
  398.         }
  399.  
  400.       MPN_POW2 (scale, 0);
  401.       assert (scalesize != 0);
  402.  
  403.       /* The number is (E - P) factors of two larger than
  404.          the fraction can represent; this is the potential error.  */
  405.       MPN_POW2 (loerr, e - p);
  406.       assert (loerrsize != 0);
  407.     }
  408.       else
  409.     {
  410.       /* The number of fractional digits is greater than the exponent.
  411.          Scale by the difference factors of two.  */
  412.       MPN_ASSIGN (r, f);
  413.       MPN_POW2 (scale, p - e);
  414.       MPN_POW2 (loerr, 0);
  415.     }
  416.       MPN_ASSIGN (hierr, loerr);
  417.  
  418.       /* Fixup.  */
  419.  
  420.       MPN_POW2 (tmp, p - 1);
  421.  
  422.       if (MPN_EQ (f, tmp))
  423.     {
  424.       /* Account for unequal gaps.  */
  425.       assert (hierrsize != 0);
  426.       cy = __mpn_lshift (hierr, hierr, hierrsize, 1);
  427.       if (cy)
  428.         hierr[hierrsize++] = cy;
  429.  
  430.       assert (rsize != 0);
  431.       cy = __mpn_lshift (r, r, rsize, 1);
  432.       if (cy)
  433.         r[rsize++] = cy;
  434.  
  435.       assert (scalesize != 0);
  436.       cy = __mpn_lshift (scale, scale, scalesize, 1);
  437.       if (cy)
  438.         scale[scalesize++] = cy;
  439.     }
  440.  
  441.       /* scale10 = ceil (scale / 10.0).  */
  442.       if (__mpn_divmod_1 (scale10, scale, scalesize, 10) != 0)
  443.     {
  444.       /* We got a remainder.  __mpn_divmod_1 has floor'ed the quotient
  445.          but we want it to be ceil'ed.  Adjust.  */
  446.       cy = __mpn_add_1 (scale10, scale10, scalesize, 1);
  447.       if (cy)
  448.         abort ();
  449.     }
  450.       scale10size = scalesize;
  451.       scale10size -= scale10[scale10size - 1] == 0;
  452.  
  453.       k = 0;
  454.       while (MPN_LT (r, scale10))
  455.     {
  456.       mp_limb cy;
  457.  
  458.       --k;
  459.  
  460.       cy = __mpn_mul_1 (r, r, rsize, 10);
  461.       if (cy != 0)
  462.         r[rsize++] = cy;
  463.  
  464.       cy = __mpn_mul_1 (loerr, loerr, loerrsize, 10);
  465.       if (cy != 0)
  466.         loerr[loerrsize++] = cy;
  467.  
  468.       cy = __mpn_mul_1 (hierr, hierr, hierrsize, 10);
  469.       if (cy != 0)
  470.         hierr[hierrsize++] = cy;
  471.     }
  472.  
  473.       do
  474.     {
  475.       mp_limb cy;
  476.       assert (rsize != 0);
  477.       cy = __mpn_lshift (r2, r, rsize, 1);
  478.       r2size = rsize;
  479.       if (cy != 0)
  480.         r2[r2size++] = cy;
  481.  
  482.       /* tmp = r2 + hierr; */
  483.       if (r2size > hierrsize)
  484.         {
  485.           cy = __mpn_add (tmp, r2, r2size, hierr, hierrsize);
  486.           tmpsize = r2size;
  487.         }
  488.       else
  489.         {
  490.           cy = __mpn_add (tmp, hierr, hierrsize, r2, r2size);
  491.           tmpsize = hierrsize;
  492.         }
  493.       if (cy != 0)
  494.         tmp[tmpsize++] = cy;
  495.  
  496.       /* while (r2 + hierr >= 2 * scale) */
  497.       assert (scalesize != 0);
  498.       cy = __mpn_lshift (scale2, scale, scalesize, 1);
  499.       scale2size = scalesize;
  500.       if (cy)
  501.         scale2[scale2size++] = cy;
  502.       while (MPN_GE (tmp, scale2))
  503.         {
  504.           cy = __mpn_mul_1 (scale, scale, scalesize, 10);
  505.           if (cy)
  506.         scale[scalesize++] = cy;
  507.           ++k;
  508.           assert (scalesize != 0);
  509.           cy = __mpn_lshift (scale2, scale, scalesize, 1);
  510.           scale2size = scalesize;
  511.           if (cy)
  512.         scale2[scale2size++] = cy;
  513.         }
  514.  
  515.       /* Perform any necessary adjustment of loerr and hierr to
  516.          take into account the formatting requirements.  */
  517.  
  518.       if (type == 'e')
  519.         cutoff += k - 1;    /* CutOffMode == "relative".  */
  520.       /* Otherwise CutOffMode == "absolute".  */
  521.  
  522.       {            /* CutOffAdjust.  */
  523.         int a = cutoff - k;
  524.         MPN_VAR (y);
  525.         MPN_ASSIGN (y, scale);
  526.  
  527.         /* There is probably a better way to do this.  */
  528.  
  529.         while (a > 0)
  530.           {
  531.         cy = __mpn_mul_1 (y, y, ysize, 10);
  532.         if (cy)
  533.           y[ysize++] = cy;
  534.         --a;
  535.           }
  536.         while (a < 0)
  537.           {
  538.         if (__mpn_divmod_1 (y, y, ysize, 10) != 0)
  539.           {
  540.             /* We got a remainder.  __mpn_divmod_1 has floor'ed the
  541.                quotient but we want it to be ceil'ed.  Adjust.  */
  542.             cy = __mpn_add_1 (y, y, ysize, 1);
  543.             if (cy)
  544.               abort ();
  545.           }
  546.         ysize -= y[ysize - 1] == 0;
  547.         ++a;
  548.           }
  549.  
  550.         if (MPN_GT (y, loerr))
  551.           MPN_ASSIGN (loerr, y);
  552.         if (MPN_GE (y, hierr))
  553.           {
  554.         MPN_ASSIGN (hierr, y);
  555.         roundup = 1;
  556.         /* Recalculate: tmp = r2 + hierr */
  557.         if (r2size > hierrsize)
  558.           {
  559.             cy = __mpn_add (tmp, r2, r2size, hierr, hierrsize);
  560.             tmpsize = r2size;
  561.           }
  562.         else
  563.           {
  564.             cy = __mpn_add (tmp, hierr, hierrsize, r2, r2size);
  565.             tmpsize = hierrsize;
  566.           }
  567.         if (cy != 0)
  568.           tmp[tmpsize++] = cy;
  569.           }
  570.       }            /* End CutOffAdjust.  */
  571.  
  572.     } while (MPN_GE (tmp, scale2));
  573.  
  574.       /* End Fixup.  */
  575.  
  576.       /* First digit.  */
  577.  
  578.       hack_digit ();
  579.  
  580.       if (type == 'e')
  581.     {
  582.       /* Exponential notation.  */
  583.  
  584.       int expt = k;        /* Base-10 exponent.  */
  585.       int expt_neg;
  586.  
  587.       expt_neg = k < 0;
  588.       if (expt_neg)
  589.         expt = - expt;
  590.  
  591.       /* Find the magnitude of the exponent.  */
  592.       j = 10;
  593.       while (j <= expt)
  594.         j *= 10;
  595.  
  596.       /* Write the first digit.  */
  597.       put (digit);
  598.  
  599.       if (low || high || k == cutoff)
  600.         {
  601.           if ((tolower (info->spec) != 'g' && prec > 0) || info->alt)
  602.         put (*decimal);
  603.         }
  604.       else
  605.         {
  606.           int stop;
  607.  
  608.           put (*decimal);
  609.  
  610.           /* Remaining digits.  */
  611.           do
  612.         {
  613.           stop = hack_digit ();
  614.           put (digit);
  615.         } while (! stop);
  616.         }
  617.  
  618.       if (tolower (info->spec) != 'g' || info->alt)
  619.         /* Pad with zeros.  */
  620.         while (--k >= cutoff)
  621.           put ('0');
  622.  
  623.       /* Write the exponent.  */
  624.       put (isupper (info->spec) ? 'E' : 'e');
  625.       put (expt_neg ? '-' : '+');
  626.       if (expt < 10)
  627.         /* Exponent always has at least two digits.  */
  628.         put ('0');
  629.       do
  630.         {
  631.           j /= 10;
  632.           put ('0' + (expt / j));
  633.           expt %= j;
  634.         }
  635.       while (j > 1);
  636.     }
  637.       else
  638.     {
  639.       /* Decimal fraction notation.  */
  640.  
  641.       if (k < 0)
  642.         {
  643.           put ('0');
  644.           if (prec > 0 || info->alt)
  645.         put (*decimal);
  646.  
  647.           /* Write leading fractional zeros.  */
  648.           j = 0;
  649.           while (--j > k)
  650.         put ('0');
  651.         }
  652.  
  653.       put (digit);
  654.       if (!low && !high && k != cutoff)
  655.         {
  656.           int stop;
  657.           do
  658.         {
  659.           stop = hack_digit ();
  660.           if (k == -1)
  661.             put (*decimal);
  662.           put (digit);
  663.         } while (! stop);
  664.         }
  665.  
  666.       while (k > 0)
  667.         {
  668.           put ('0');
  669.           --k;
  670.         }
  671.       if ((type != 'g' && prec > 0) || info->alt)
  672.         {
  673.           if (k == 0)
  674.         put (*decimal);
  675.           while (k-- > -prec)
  676.         put ('0');
  677.         }
  678.     }
  679.     }
  680.  
  681. #undef    put
  682.  
  683.   /* The number is all converted in BUF.
  684.      Now write it with sign and appropriate padding.  */
  685.  
  686.   if (is_neg || info->showsign || info->space)
  687.     --width;
  688.  
  689.   width -= bp - buf;
  690.  
  691.   if (!info->left && info->pad == ' ')
  692.     /* Pad with spaces on the left.  */
  693.     while (width-- > 0)
  694.       outchar (' ');
  695.  
  696.   /* Write the sign.  */
  697.   if (is_neg)
  698.     outchar ('-');
  699.   else if (info->showsign)
  700.     outchar ('+');
  701.   else if (info->space)
  702.     outchar (' ');
  703.  
  704.   if (!info->left && info->pad == '0')
  705.     /* Pad with zeros on the left.  */
  706.     while (width-- > 0)
  707.       outchar ('0');
  708.  
  709.   if (fwrite (buf, bp - buf, 1, s) != 1)
  710.     return -1;
  711.   done += bp - buf;
  712.  
  713.   if (info->left)
  714.     /* Pad with spaces on the right.  */
  715.     while (width-- > 0)
  716.       outchar (' ');
  717.  
  718.   return done;
  719. }
  720.  
  721. #ifndef NDEBUG
  722. static void
  723. mpn_dump (str, p, size)
  724.      const char *str;
  725.      mp_limb *p;
  726.      mp_size_t size;
  727. {
  728.   fprintf (stderr, "%s = ", str);
  729.   while (size != 0)
  730.     {
  731.       size--;
  732.       fprintf (stderr, "%08lX", p[size]);
  733.     }
  734.   fprintf (stderr, "\n");
  735. }
  736. #endif
  737.